拉普拉斯矩阵与图卷积神经网络

ref1. https://zhuanlan.zhihu.com/p/84649941

ref2. https://zhuanlan.zhihu.com/p/84649941

对于实数域，导数 \nabla f = \sum_{i=1}^K \frac{df}{dx_i}\vec{e_i} ，散度 div f = \sum_{i=1}^K \frac{df}{dx_i} .

拉普拉斯算子 \Delta f=\nabla^2 f=\sum_{i=1}^K \frac{d^2f}{dx_i^2} ，含义为向量场中的梯度净流出程度，即导数的散度；

考虑扩展到图领域。

定义 G=(V,E) ， |V|=N, e_{i,j}\in E ， [e_{i,j}\in E]\triangleq (i\to j) ， |E|=M 。

为了方便表示，设 e_{i,j}=\{0,1\} ，其中1表示 e_{i,j}\in E ，反之则不存在。令边权为 w_{i,j} 。

设每个节点为 x_i ，每个节点的特征为 f(x_i)\in R (高维同理)， \vec{f}=(f(x_1)...f(x_N))^T

.

首先定义“边的梯度”：

令 T\in R^{N\times M} 为权值矩阵，若 e_{i,j}\in E ，且 e_{i,j} 对应第 t 条边，则 T_{i,t}=w_{i,j},T_{j,t}=-w_{i,j} .

设 e_{i,j}=1 表示“正交”，反之“不正交”。

则类似于实数域 \frac{df}{dx}=\frac{f(x+dx)-f(x)}{dx} ，将 dx \triangleq w_{i,j} ，有 \nabla w_{i,j} = T_{i,j}(f(x_j)-f(x_i)) .

令边的梯度向量为 Y\in R^{M} ，则 Y=T^T \vec{f}

.

接着定义“点的拉普拉斯”：

类似于实数域，拉普拉斯即导数的散度，“点的散度”即“边梯度之散度”。

令 H\in R^{N\times M} 为连接矩阵，若 e_{i,j}\in E ，且 e_{i,j} 对应第 t 条边，则 H_{i,t}=1,H_{j,t}=-1 .

则”以点为对象”的”边梯度的散度”等于将 dx \triangleq H_{i,j} 后，对”边梯度”求梯度后再求和。

\nabla (\nabla W_{i,j}) = H_{i,j} \nabla w_{i,j} 。故对于节点而言， \Delta f(x_i) = \sum_{v\in N(i)} H_{i,j}\nabla w_{i,j} 。

令节点的拉普拉斯向量为 Z\in R^{N} ，则 Z=HY=HT^T\vec{f}

.

其中 L=HT^T 即为拉普拉斯矩阵。当 w_{i,j}=e_{i,j} ，则有 L=HH^T ，其中 H\in R^{N\times M} .

不难发现 H_{u,v}W_{u,v}=w_{u,v} 恒成立，故拉普拉斯矩阵与边的方向无关。

因此，拉普拉斯矩阵一般用于描述无向图。

令 W\in R^{N\times N}, A\in R^{N\times N} ，其中若 e_{i,j}\in E ，则 A_{i,j}=A_{j,i}=1, W_{i,j}=W_{j,i}=w_{i,j} .

我们称 L=WA^T=HT^T 为拉普拉斯矩阵。

.

下面抽离出来看下拉普拉斯到底求解了什么： Z_u = \sum_{v\in N(u)} H_{u,v}W_{u,v}(f(x_u)-f(x_v)) .

根据对无向图性质的观察，则 Z_u = \sum_{v\in N(u)}w_{u,v}(f(x_u)-f(x_v)) 。

故 Z_u = f(x_u)\sum_{v\in N(u)} w_{u,v} - \sum_{v\in N(u)} w_{u,v}f(x_v) 。

令 d_u = \sum_{v\in N(u)} w_{u,v}，D=Diag(d_1...d_N) ，则 Z = D\vec{f} - W \vec{f}=(D-W)\vec{f} .

因此，拉普拉斯矩阵有： L=D-W 。当 G 为无权图时， L=D-A 。

更多性质见：ref3. https://mathweb.ucsd.edu/~fan/research/cb/ch1.pdf

谱分解即特征分解。由于 L 为半正定矩阵，故 L 一定可分解： L=U\Lambda U^T .

类似的， A 和 W 也为对称矩阵，故一定可以谱分解.

我们直到，GNN的一般形式为： f^{(k+1)}(x_u)=\sigma(f(x_u)^{(k)} + MERGE_{v\in N(u)}(f(x_v)^{(k)})) .

使用均值融合并用 A 传播，则： \vec{f}^{(k+1)} = \sigma(A\vec{f}^{(k)}F) ，其中 A 用于传播， F 用于transform.

可以看到， A^k 即在图上进行 k 次传播。

而令 A=U\Lambda U^T ，则 A^k = U \Lambda^k U^T ，而 \Lambda^k 的计算即使用蛮力法也是 O(Nk) 的。

A^k \to \Lambda^k 即将“时域的卷积”转化为”频域的乘法”。”频域”的含义见下文。

又称 U^T 为傅里叶变换，即 \hat{f}=U^Tf 。则傅里叶逆变换为 f=U\hat{f} 。

类似于主成分分析(PCA)，谱分解后可以只取若干部分， \lambda 值越小意味着该成分”信息越少”。

ref4. https://zhuanlan.zhihu.com/p/551528542

ref5. https://zhuanlan.zhihu.com/p/600232548

拉普拉斯矩阵 L 的谱分解具有特殊的性质，其谱分解的特征值 \lambda 是具有实际意义的”频域”，其含义为图的局部波动性（注意 \lambda_i \geq 0 ，对应波动性大于零）。简单来说，由于 u_i \in R^N ，把 u_{i,j} 分配到$j$节点作为权重，则 \lambda_i 越大，局部的波动性越大。可视化结果见ref4.

首先对 L 的半正定性进行证明，半正定性只需要证明二次型大于等于零。

f^TLf = f^T(D-W)f = \sum_{i}d_i f_i^2 - \sum_{i,j} w_{i,j}f_if_j = \sum_{i}(d_i-w_{i,i})f_i^2 - \sum_{i,j}w_{i,j}f_if_j .

重新分配后， f^TLf = \sum_{e_{i,j}\in E}w_{i,j}(f_i-f_j)^2 \geq 0 .故 L 半正定。

一个有意思的事情是，考虑无权图 G ，定义特征 f 的波动率为 VT(f)=\sum_{e_{i,j}\in E}(f_i-f_j)^2 .

则根据上述证明， VT(f)=f^TLf 。

令 f=u_i ，则 VT(u_i)=u_i^T L u_i = u_i^T (U\Lambda U^T)u_i = (U^Tu_i)^T \Lambda (U^T u_i)=\lambda_i .

故 \lambda_i = \sum_{e_{x,y}\in E}(u_{i,x}-u_{i,y})^2 ，即对应以 u_i 为权重的波动率。

.

另外，由于 \lambda_1=0, u_1=(1,1...1)^T ，故 u_1 对应完全没有波动，每个节点权重都相等。

由于 u_i^Tu_1=I(i=1) ，故对于 i\neq1 ，有 \sum_{j=1}^N u_{i,j}=0 .

即对于 i\neq 1 ， MEAN(u_{i})=0 ， VAR(u_{i})=\lambda_i 。

实证上，通过仿真可以看到， u_i 穿越零点的次数和 \lambda_i 基本成正比。

.

因此， u_i 和 \lambda_i 某种程度上对应滤波，即较小的 \lambda_i 对应的 u_i 关注低通部分，即图的全局信息，而较大的 \lambda_i 对应的 u_i 关注高通部分，即图的局部信息。”频域”在图中对应局部波动性。

ref6. https://blog.csdn.net/yyl424525/article/details/100058264

ref7. Simplifying Graph Convolutional Networks**.** Felix Wu, Tianyi Zhang, Amauri Holanda de Souza Jr., Christopher Fifty, Tao Yu, Kilian Q. Weinberger. ICML 2019

为了数值稳定，存在一些常见的标准化版本的拉普拉斯矩阵。

下面将情况进一步泛化，设 \vec{f} 为特征，则基于拉普拉斯的消息传递为： L\vec{f}=U\Lambda U^T \vec{f}

如上文所述， U^T 是傅里叶变换矩阵， u_i 本身包含了”频域信息”， \Lambda_i 只是显式的对波动进行了显示。

因此，实际上 \Lambda 可以替换为其它”卷积核”，设为 (g_1,g_2...g_N) ，令 G=Diag(g_1,...g_N) .

则 \hat{f} * \hat{g} = (U^T g)\circ (U^T f) ，其中 \circ 为Hadamard，即元素点乘。

故一次卷积操作对应： f*g = U(\hat{f}*\hat{g})=U((U^T f)\circ (U^Tg))=UGU^T f .

基于Spectral的神经网络主要是操作 g=(g_1,...g_N) 。

.

一个例子是ChebNet： G=\sum_{k=0}^K \theta_k \Lambda^k 。其中 \Lambda 是 L_{sys} 的谱分解频率。

对ChebNet进行一阶近似即GCN。

令 \theta_0=1, \theta_1=-1, \theta_{k>1}=0 ，则 GCN=f*g=(I+D^{-\frac{1}{2}}AD^{-\frac{1}{2}})f .

但GCN的作者使用了Renormalization的技巧，令 \tilde{A} = I+A ， \tilde{D} = I+D ，

则 RGCN=\tilde{D}^{-\frac{1}{2}} \tilde{A} \tilde{D}^{-\frac{1}{2}}f .

下面考虑对 GCN：I+D^{-\frac{1}{2}}AD^{-\frac{1}{2}} ， RGCN：\tilde{D}^{-\frac{1}{2}} \tilde{A} \tilde{D}^{-\frac{1}{2}} 进行谱分析。

设 L = (D-A) = U \Lambda U^T ， \tilde{L} = (\tilde{D}-\tilde{A}) = \tilde{U} \tilde{\Lambda} \tilde{U}^T ，

则 GCN: g_i^K = (2-\lambda_i)^K ， \tilde{g_i}^K=(1-\tilde{\lambda_i})^K .

由于 \lambda_i \in [0,2] ，故 g_i 是低通滤波器，但当 \lambda_i < 1 时， g_i 存在放大作用，这可能导致过平滑。

为了更好的分析 \tilde{g_i} ，我们先对 D^{-\frac{1}{2}}A D^{-\frac{1}{2}} 进行谱分析，其对应的 g_i^{norm} = (1-\lambda_i)^K 。

对于 g_i^{norm} 滤波器，其不存在放大过平滑的问题，但当 \lambda_i >1 时且 K 奇数时，存在负数核的情况。

但存在定理(见ref7)： 0=\lambda_1 = \tilde{\lambda_1}